Spatial correlations of trapped ID bosons in an optical lattice 
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We investigate a quasi-one dimensional system of trapped cold bosonic atoms in an optical lattice 
by using the density matrix renormalization group to study the Bose-Hubbard model at T = for 
experimentally realistic numbers of lattice sites. It is shown that a properly rescaled one-particle 
density matrix characterizes superfiuid versus insulating states just as in the homogeneous system. 
For typical parabolic traps we also confirm the widely used local density approach for describing 
correlations in the limit of weak interaction. Finally, we note that the superfiuid to Mott-insulating 
transition is seen most directly in the half width of the interference peak. 
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During the last years enormous progress was made in 
the experimental manipulation of cold atoms in optical 
lattices. Recently, Greiner et al. succeeded in driv- 
ing a transition between a superfiuid (SF) and a Mott- 
insulating (MI) state in a system of ultracold bosonic 
atoms in an optical lattice as predicted by Jaksch et al. 
0. In contrast to solid state realizations the experi- 
mental setup involves the application of an additional 
parabolic trapping potential that causes a state in which 
the two phases, though spatially separated, coexist 0. 
Due to the inhomogeneity the usual characterization of 
the SF to MI transition by the asymptotic behaviour of 
the one-particle density matrix does not apply. Moti- 
vated by this, we use the density matrix renormalization 
group (DMRG) H to study how the parabolic confining 
potential influences the one-particle density matrix and 
its Fourier transform, which is related to the interference 
pattern observed in the experiments || . We find that by 
a simple rescaling, the decay of the correlations can be 
used to characterize the occuring states, just as in the 
homogeneous case. We further confirm the applicabil- 
ity of the standard local density approximation to the 
inhomogeneous system fj| for weak interactions by com- 
paring it to the DMRG results for the correlation func- 
tions. Studying experimentally accessible quantities we 
find that the half width of the interference peak contains 
the essential information about the state of the system. 

Model: Ultracold bosonic atoms in an optical lattice 
@ can be described by a Bose-Hubbard model 
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where b\ and bj are the creation and annihilation oper- 

3 J 

ators on site j and hj — bjbj is the number operator 
7J. This Hamiltonian describes the interplay between 
the kinetic energy due to the next-neighbour hopping 
with amplitude J and the repulsive onsitc interaction U 
of the atoms. By tuning the lattice depth in the exper- 
iment, the parameter u = U j J can be varied over sev- 
eral orders of magnitude. To investigate the properties 
of the ID Bose-Hubbard model, we apply the DMRG, a 
quasi-exact numerical method, very well suited to study 



strongly correlated quasi ID quantum systems with a 
large number of sites at zero temperature It has 
been successfully applied to spin, fermionic and bosonic 
quantum systems including the homogeneous [jj and the 
disordered [ljj Bose-Hubbard model. We used the fi- 
nite size DMRG algorithm [|| which is better suited for 
an inhomogeneous system, since it gives the system the 
possibility to evolve further after the final length of the 
system is reached. Additionally some tricks are applied 
to circumvent problems which arise due to the sparse fill- 
ing at the boundaries. The numerical results were tested 
to be convergent in the cut-offs used for the length of the 
system, the number of states kept for the Hilbert space, 
and the number of states allowed per site. Uncertainties 
given below are determined by comparing data of differ- 
ent parameter sets. 

State diagram: The confining trap of the experiment 
which consists of a magnetic trap and the confining 
component of the laser which generates the optical lat- 
tice, can be modeled by setting £j — V£ ap (a(j — jo)) 2 
in Eq. (JJJ, where a is the lattice constant. We choose 
the strength of the trap proportional to the onsite in- 
teraction, i.e. V^r ap = vqU, since this guarantees that 
when the optical lattice depth, corresponding to the pa- 
rameter u in the Bose-Hubbard model, is changed, the 
size of the system does not vary much for a fixed particle 
number. This is consistent with the experimental real- 
ization, in which the total size of the condensate is essen- 
tially independent of the lattice depth. In the presence 
of a parabolic trap at average filling of approximately 
one particle per site, one can distinguish three states of 
the system (see |2,|s|): (a) for u < u c i, the particle oc- 
cupancy is incommensurate over the whole system; (b) 
for u c i < u < u C 2, regions with incommensurate and 
commensurate occupancy coexist; and (c) for u > u C 2, 
the main part of the system is locked to commensurate 
filling and only at the boundaries small incommensurate 
regions exist. For small particle numbers, state (b) does 
not occur. A sketch of the state diagram is presented in 
Fig. ^ (A) • The insets show the characteristic shape of 
the particle distribution for the three states. For state 
(b) the exact locations of the interface between the com- 
mensurate and the incommensurate regions are difficult 
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FIG. 1: A) Sketch of the state diagram for v = 4/64 2 . 
The insets sketch the shape of the density distribution in the 
states. (B) Sketch of the phase diagram of the homogeneous 
system: chemical potential p, versus 1/u. The different sym- 
bols in (B) mark the locations of the chemical potential values 
in the local density approximation that correspond to the lo- 
cations in the density profiles marked in (A). 



to determine. This is due to the fact that these sites cor- 
respond in the homogeneous system to the critical pa- 
rameter regimes at the phase transition, where strong 
fluctuations and extreme sensitivity to boundary condi- 
tions make a numerical investigation very difficult. 

Reseated correlations: To get a better understanding 
of the three states (a)-(c), we study the properties of the 
rescaled one-particle density matrix, 



(2) 



in which the leading density dependence of bj oc */nJ is 



divided out. In the absence of density fluctuations Cj(r) 
is just the pure phase correlation function (e i ^ j e~ i ^ j+r ). 
At the two particle level, the equivalent step is going 
from the two-particle density p( 2 \x\, x^) to the dimcn- 
sionless two particle distribution function g^> (aq, X2) — 
p^(xi,X2)/ ' p^{xi)p^{x2). Remarkably, we find that 
by this simple rescaling, the signatures of the SF and MI 
phases in the homogeneous system, namely an algebraic 
or exponential decay, Cj(r) oc A\r\~ K / 2 and cx Be^^^, 
respectively, can be recovered approximately even in the 
presence of a parabolic confining potential. For weak 
interactions, u < u c \, [Fig. (a)] Cj(r) decays approx- 
imately algebraically with r.In the intermediate regime, 
u c i < u < u c2 , [Fig. [2(b)] the decay in the regions where 
the density is incommensurate is still algebraic, whereas 
in the regions where the density is locked, it shows an ex- 
ponential behaviour. Increasing the interaction further, 
u > u C 2, [Fig. [3(c)] the incommensurate regions disap- 
pear and the correlations decay exponentially. 

Hydrodynamical approach: It is instructive to com- 
pare the numerically exact DMRG results to a hydrody- 
namical treatment of the interacting ID Bose gas [llj 
combined with a local density approximation. In the 
hydrodynamical approach the low-energy fluctuations of 
the system are described by two conjugate fields, the 
phase fluctuations (j)(x) and the density fluctuations 9{x). 
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FIG. 2: Scaled correlations Cj(r) [Eq. [5] for different fixed 
sites j are plotted as a function of r for different values of u. 
For the coexistence region (b) a shallower trapping potential 
is chosen, such that the extents of both the incommensurate 
and the commensurate region are large enough to allow iden- 
tification of the algebraic and exponential behaviour. 



This approach can be generalized to the case of inhomo- 
geneous systems |(J by taking the density fluctuations 
around a smooth, spatially dependent density profile 
n(x). An equivalent procedure was used for ID Fermionic 
gases by Recati et al. [l^. The Hamiltonian becomes 



H = — J dx {v j (x)(d x (j)) 2 + v N (x)[d x 6 - Trn(x)) 2 } , 

precisely as in the homogeneous case, except that 
n(x), and therefore Vj(x) = Trhn(x)/m and v^(x) — 
(7T?i) _1 (lj^)| n=n ( x ), now depend on x. To account for the 
inhomogeneity, the local density approximation p[n(x)] + 
V(x) = /i[n(0)] was used to obtain the mean density 
profile [l3j. Based on this approximation Gangardt and 
Shlyapnikov have shown that the normalized matrix 
elements of the one-particle density matrix are given by: 



C(x) 



\2x\ 



y/n{x)n{-x) \h{x) 
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(3) 



where K is the exponent and l c the longitudinal correla- 
tion length. Eq. 10 is derived assuming \2x\ l c . Spe- 
cializing to weak interaction, i.e. 7 = 1/dn <C 1, the ap- 
proximations l c (x) w \J~~^e~) an< ^ K{x) w 1/ (tt-\J dn(x)^j 
hold, where d cx l\/azD is the characteristic length of 
the interaction, d depends on the 3D scattering length 
a^o and the amplitude l± of the transverse zero point 
oscillation. The condition \2x\ ^> l c breaks down at the 
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FIG. 3: Quasi-exact DMRG results for C(j) (symbols) are 
compared to equation obtained by the hydrodynamical 
approach |(| (lines). We used n(x) = no(l — (x/R) 2 ), where 
no and R are determined by fitting to the DMRG results 
(see insets). The uncertainties are obtained by varying the fit 
range in the sensible region away from the boundaries. 



boundaries, where n(x) vanishes causing a divergence in 
l c (x). Comparing [Eq. to the quasi-exact results of 
DMRG, it turns out that the local density approach de- 
scribes very well the rescaled correlations in the inhomo- 
geneous systems for 7 < 2. To this end we fitted the 
function C(x) [Eq. ©] to the corresponding DMRG re- 
sults, using only d as fitting parameter [Fig. [3] ■ We find 
very good agreement in the bulk of the SF regions in 
both, the purely SF state [Fig.|3|(a)] and the coexistence 
state [Fig. El (b)]. The quality of the agreement is some- 
what surprising, because the pure state (7 = 0.6) and 
the coexistence state (7 = 1.7) are in an intermediate 
regime between the Thomas- Fermi limit (7^1) and the 
Tonks gas (7 3> 1), where the density profile is no longer 
parabolic |13| . 

Interference pattern: We investigate how the informa- 
tion contained in the interference pattern is influenced 
by the confining potential. If the interaction between the 
atoms after switching off the confining potentials is weak, 
i.e. E pot <C £"kin, the measured absorption images reflect 
the momentum distribution obtained from the Fourier 
transform of the one-particle density matrix [{| 
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where M is the number of sites in the chain and N the to- 
tal number of particles. For the parameters studied here, 
the approximation of a negligible contribution of the in- 
teraction energy to the time of flight images is valid for 
all momenta in the second or in higher Brillouin-zones. 
Indeed, these momenta are of order 2hirs/L, where s € N 

and s > M. Thus f=t cx T¥t^t¥^ °c " « KT 2 

Skin (7Tfts / LV J (2m) a 

for n^D < 1.5/a 3 and a s /a like in The function p(k ) 
has been studied for very small systems numerically [14j , 
with the hydrodynamical approach 0] for a ID homo- 
geneous system and for the confined system in 3D [j| and 
ID [l3|. In Fig. H we plot the DMRG results (symbols) 
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FIG. 4: Interference pattern for the system with (A) open 
boundaries and with (B) parabolic trap for different values 
of u. Symbols are the results of the DMRG (maximal un- 
certainty 0.1) and lines the results of the approximations ex- 
plained in the text. The insets enlarge the scale of the y-axis. 
For a homogeneous system u c (n = 1) ~ 3.37 is the critical 
value in the thermodynamic limit according to [ij. 



for the function p(k) for several values of the parame- 
ter u, comparing the homogeneous system (e, = 0) with 
open boundary conditions (A) to the parabolic system 
(B). In the homogeneous system with commensurate fill- 
ing, n = N/M—l [Fig. 0] (A)] we find a very sharp peak 
at small momenta for u < u c . If u is increased, the peak 
height decreases smoothly. The half width w [Fig. [S] (A)] , 
however, shows a clear upturn. This upturn signifies a 
phase transition, since it stems from the behaviour of 
the correlation length £ (oc u> -1 ), which diverges in the 
SF phase (£ cx L) and becomes finite in the MI phase 
(£ cx A -1 , where A is the energy gap). For the parabolic 
system [Fig. 0] (B)], the interference pattern for small 
and large u is similar to the interference pattern in the 
homogeneous system. In the intermediate regime, how- 
ever, it shows a more complex behaviour, which is most 
clearly evident in w [Fig. [S] (B)]. For small particle num- 
bers (JV = 40), w is very small for u < u c \ and rises 
continuously for u > u c \. In contrast, for larger particle 
numbers (JV = 50, 60) three different regimes correspond- 
ing to the three different states in Fig.QJare observed: (a) 
for u < Ud, w is very small, (b) for u c i < u < u C 2, w 
rises slowly, until at u ~ u C 2 it shows a sudden jump-like 
increase, (c) for u > u C 2, it continues to rise strongly. 
That means that in the SF (a) and the MI (c) state the 
behaviour of w resembles that of the homogeneous sys- 
tem. This is as expected, since the rescaled correlations 
show the same decay as in the corresponding homoge- 
neous phases. In the intermediate regime (b), however, 
it shows a new behaviour, a slow increase, which is due 
to the coexistence of the SF and the MI state. The SF 
region determines mainly the height of the interference 
peak, while its broadening is due to the presence of the 
MI region. In the crossover region between the totally 
incommensurate and the coexistence region, the inter- 
ference pattern shows additional ocillations with period 
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FIG. 5: Half width of the interference peak for the homo- 
geneous (A) and the parabolic (B) system. Arrows in (A) 
mark the critical value of u c in the thermodynamic limit (solid 
and dashed for n = 1 and n = 2 respectively) according 
to Arrows in (B) mark the three different regimes de- 
scribed in the text. To relate u to the corresponding lattice 
depth VL — Viax/Er of experiments, we assumed that the 
depths in the two perpendicular dimensions where fixed to 
Vut,x/E r = 50. 

2ir/l, where I is the distance between the two outer SF 
regions, due to the appearance of relatively strong corre- 
lations between the latter. Similar oscillations were seen 
in [Hi ] . In smaller systems such as in the effect is more 
pronounced causing well-separated satellite peaks . 

Finally, let us investigate to what extent the properties 
of the interference patterns in Fig. 0] can be understood 
in terms of simple phenomenological approximations for 
{b'jbj,} in the homogeneous and the rescaled correlations 
Cj (r) in the inhomogeneous system. Once the character- 
istic quantities K and £ have been identified (in this case 
by fitting to DMRG results), our simple rescaling pro- 
cedure captures most of the essential observable physics. 
To illustrate this we show in Fig. 0] (A) in addition to the 
DMRG results, results (lines) obtained by approximat- 
ing (b]bj,) in Eq. Q by A\j - j'\- K/2 and Be'^'l/^ 
for small and large u, respectively. The values of K and 
£ are determined by fitting (b^b^,) to DMRG results (not 
shown here) . The constants A and B are chosen such that 
the value at k — agrees with the DMRG results. In Fig. 
51(B) the approximation (lines) are obtained analogously 
by taking the density scaling into account, i.e. replacing 
(bjbj,) by the algebraically and the exponentially decay- 
ing functions times the scaling factor ^JnJrTf. We use the 
density distribution rij = no(l — (j — jo) 2 /R 2 ) for u = 1, 
and rij = 1 for u = 9. The parameters K and £ are 
determined by fitting the rescaled correlation functions. 
Comparing the DMRG data to the approximation we see 
in Fig. 0] that this simple approximation works very well 
for small values of feet; in particular, it reproduces the 
correct shape of the peak [even including the small non- 
monotonities which are due to the finite sum in Eq. (0}]. 
This underlines that p{k) is mainly determined by the 
decay of the (un)scaled correlations. Clearly our calcula- 



tions in ID cannot be compared quantitatively with the 
experiments in a 3D lattice 0. Recently, however, an ar- 
ray of truly ID Bose systems has been created ^jj- With 
an additional lattice potential our predictions can then 
be tested quantitatively |l8| . In the experimental real- 
ization one typically has several ID systems next to each 
other with different particle number, hence the location 
of the sharp upturn in the half width [Fig. E1(B)], will be 
smeared out, since the critical value u C 2 depends on the 
particle number. Nevertheless, we expect in particular 
the strong, jump-like increase between the coexistence 
state and the MI state to remain observable. 

In conclusion, we have found that the correlation func- 
tions of a parabolically confined system, after a remark- 
ably simple rescaling, show approximately the familiar 
algebraic and exponential behaviour of the SF and MI 
phases in the homogeneous system. We investigated as 
well the applicability of the local density approximation 
in a parabolic system in the limit of weak interaction 
and find good agreement with the DMRG results. More- 
over, if the experimental system consists of ID tubes with 
almost the same average filling, the half width of the in- 
terference peak can be used to distinguish the different 
types of states that occur experimentally. 
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